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ABSTRACT 

We study the evolution of an isolated, spherical halo of self-interacting dark matter (SIDM) in the gravothermal 
fluid formalism. We show that the thermal relaxation time, t^, of a SIDM halo with a central density and velocity 
dispersion of a typical dwarf galaxy is significantly shorter than its age. We find a self-similar solution for the 
evolution of a SIDM halo in the limit where the mean free path between collisions. A, is everywhere longer than 
the gravitational scale height, H. Typical halos formed in this long mean free path regime relax to a quasistationary 
gravothermal density profile characterized by a nearly homogeneous core and a power-law halo where p oc r"^ 
We solve the more general time-dependent problem and show that the contracting core evolves to sufficiently high 
density that A inevitably becomes smaller than H in the innermost region. The core undergoes secular collapse 
to a singular state (the "gravothermal catastrophe") in a time tcou ~ 290f,., which is longer than the Hubble time 
for a typical dark matter-dominated galaxy core at the present epoch. Our model calculations are consistent with 
previous, more detailed, N-body simulations for SIDM, providing a simple physical interpretation of their results 
and extending them to higher spatial resolution and longer evolution times. At late times, mass loss from the 
contracting, dense inner core to the ambient halo is significantly moderated, so that the final mass of the inner core 
may be appreciable when it becomes relativistic and radially unstable to dynamical collapse to a black hole. 

Subject headings: cosmology:theory — dark matter — galaxies:formation — dynamics 



1. INTRODUCTION 

Cold dark matter has become the paradigm for explaining ob- 
servations of large scale structure of the universe. Flat cosmo- 
logical models composed of cold dark matter (CDM) and a cos- 
mological constant (or quintessence) combined with a nearly 
scale-invariant, adiabatic spectrum of density fluctuations pro- 
duce excellent fits to observed structure on large (^ 1 Mpc) 
scales (Bahcall et al. 1999). This class of models - collectively 
labeled as ACDM - has received substantial additional sup- 
port from the recent cosmological observations of supernovae 
at high redshift (Riess et al. 1998; Perlmutter et al. 1999) and 
the fluctuations in the cosmic microwave background (Netter- 
field et al. 2001; Halverson et al. 2001). The nature and proper- 
ties of the dark matter remain unknown, however, and present 
one of the greatest challenges of current cosmological research. 

A recent suggestion by Spergel & Steinhardt (2000) that ob- 
servations on galactic and subgalactic scales offer unique clues 
regarding the properties of the dark matter has attracted much 
attention. The basic argument is that the "standard" CDM mod- 
els, where dark matter particles are assumed to interact only 
through gravity, may be in conflict with some observed fea- 
tures on smaller (< 1 Mpc) scales. In particular, CDM models 
(Navarro et al. 1997; Moore et al. 1999) tend to predict density 
profiles of dark matter halos which are cuspy and have cen- 
tral densities of 1 Mq pc"^ and larger, whereas observations 
of dwarf and low surface brightness galaxies seem to favor 
the presence of a relatively flat core, with typical densities of 
0.02 Mq pc"^ (Firmani et al. 2001). Other apparent discrepan- 



cies from the simulations of CDM are that they tend to predict: 
a) more than ten times as many dwarf galaxies as are observed 
for the Local Group, b) disks which are small and have too lit- 
tle angular momentum, c) triaxial clusters instead of relatively 
spherical ones; see Dave et al. (2001) for discussion and ref- 
erences. The thrust of the suggestion by Spergel & Steinhardt 
(2000) is that these possible discrepancies may be alleviated if 
dark matter particles have some additional type of interaction, 
besides gravity, that affects structure formation only on these 
smaller scales^. 

The physical requirement for such "self-interacting dark mat- 
ter" (SIDM) is that on the smaller scales, the typical density of 
dark matter particles is large enough so that interactions have 
a nonnegligible effect, while on large scales the probability of 
interactions is low enough so that the favorable properties of 
CDM vis a vis observed structure are maintained. Spergel & 
Steinhardt (2000) speculated that if the dark matter particles 
scatter off one another through this interaction, the entropy of 
their phase space must increase while their trajectories will be- 
come isotropic rather than radial, processes which will lead to a 
dark matter core with a shallower density profile. The ansatz is 
that a dark matter particle in a typical galactic halo should un- 
dergo several collisions per Hubble time to achieve such halo 
structure. Since the relevant mean free path is 1 - 1000 kpc 
and typical dark matter densities (e.g, in the solar neighbor- 
hood) are < 1 GeV cm"-^, the cross section for SIDM is required 
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(2001) concluded from numerical simulations that the preferred 
range for the cross section is about 0.5 - 5 cm^ gm~' . Inter- 
estingly, this coincides with the range of typical low energy 
hadronic physics, implying that a natural candidate for a SIDM 
particle could be a light, supersymmetric hadron (Wandelt et al. 
2000; Balberg, Farrar & Piran 2001). 

Following the suggestion by Spergel & Steinhardt (2000), 
several numerical studies of SIDM were conducted via N-body 
simulations, amended to include a Monte-Carlo algorithms for 
particle-particle scattering. Yoshida et al. (2000) and Dave 
et al. (2001) examined SIDM structure formation in a universe 
with an initial density fluctuation spectrum. Burkert (2000) and 
Kochanek & White (2000) studied the evolution on a single, 
isolated halo, with initial conditions based on the Hernquist 
(1990) model, which describes a noninteracting CDM halo (this 
choice is probably inappropriate for SIDM - see below, and also 
in Dave et al. 2001). While there is significant disagreement 
about whether various aspects of the results are favorable for 
the SIDM model, the general trend in the N-body simulations 
is that the centers of halos do indeed become flatter with lower 
densities with respect to standard ACDM models. 

In this work we revisit the issue of the density profile of 
a SIDM halo by studying the evolution of such a halo via 
a gravo thermal fluid approximation. In the gravo thermal ap- 
proach, the ensemble of gravitating particles is approximated 
as a fluid in quasi- static virial equiUbrium. The effective tem- 
perature is identified with the square of the one dimensional ve- 
locity dispersion, and thermal heat conduction is employed to 
reflect the manner in which orbital motion and scattering com- 
bine to transfer energy in the system. This formaUsm was origi- 
nally introduced for the study of globular star clusters (see, e.g, 
reviews by Lightman & Shapiro 1978 and Spitzer 1987), and 
has proven to be very successful in understanding the evolution 
of these systems (Larson 1970; Hachisu et al. 1978; Lynden- 
Bell & Eggleton 1980). This agreement comes about despite 
the fact that star clusters are only weakly collisional and have 
long collision mean free paths greatly exceeding the size of the 
cluster, and thermalization is achieved by the cumulative ef- 
fect of repeated distant, smaU-angle, gravitational (Coulomb) 
encounters. In fact, the gravothermal fluid description is much 
better suited to SIDM halos where the thermalizing particle in- 
teractions are close encounter, large-angle (hard sphere) scat- 
terings. It is reassuring, nevertheless, that even in the case 
of weakly collisional systems like star clusters, the gravother- 
mal fluid model does reproduce many of the results found from 
more fundamental analyses of the collisional Boltzmann equa- 
tion, with colhsions treated in the Fokker-Planck approxima- 
tion (Henon 1961; Spitzer and Hart 1971; Cohn 1979, 1980; 
Marchant & Shapiro 1980). 

An isolated self-gravitating cluster of particles in virial equi- 
librium will relax via collisions to a state consisting of an ex- 
tended halo with a shallow temperature gradient, surround- 
ing a hotter, central core region. As time advances, the core 
transfers mass and energy through the flow of particles and 
heat to the extended halo. The thermal evolution timescale 
of the dense core is much shorter than that of the extended 
halo, which essentiaUy serves as a static heat sink. As the 
core evolves it shrinks in size and mass, while its density and 
temperature grow. Increase of central temperature induces fur- 
ther heat transfer to the halo, leading to a secular instabil- 
ity on a thermal (colUsional relaxation) timescale. The secu- 
lar evolution of the core towards infinite density and tempera- 



ture but zero mass is known as "gravothermal collapse" or the 
"gravothermal catastrophe" (Lynden-Bell & Wood 1968; see 
reviews in Lightman and Shapiro 1978 and Spitzer 1987). A 
critical modification to the gravothermal scenario is that a dy- 
namical instability occurs when the particle velocities in the 
core, or, equivalently, when the central potential, becomes rel- 
ativistic. It is well known that fluid stars (i.e. highly collisional 
gases) in hydrostatic equilibrium are unstable to radial col- 
lapse on a dynamical timescale when they are sufficiently com- 
pact and relativistic. Zel'dovich & Podurets (1965) originally 
conjectured and Shapiro and Teukolsky (1985a,b, 1986) sub- 
sequently demonstrated that colUsionless gases in virial equi- 
librium also experience a radial instability to collapse on dy- 
namical timescales when their cores are sufficiently relativistic. 
In either case, this dynamically instability will terminate the 
epoch of secular gravothermal contraction and will lead to the 
catastrophic collapse of a finite mass core, which must end as a 
black hole. Simulations of the catastrophic collapse of relativis- 
tic coUisionless clusters have been performed recently, in part 
to explore the possibility of creating the supermassive black 
holes (SMBHs) inferred to exist in most galaxies and quasars 
(Shapiro & Teukolsky 1985c); see Shapiro & Teukolsky (1992) 
for a review. 

In the case of globular star clusters and gravitational scatter- 
ing, multiple, distant, small-angle scattering events dominate 
over the very rare, close encounter, large-angle scatterings until 
the core is very dense (Lightman & Shapiro 1978). Goodman 
(1983) showed that the last stage of secular core collapse, when 
large-angle scatterings in the core are not neghgible, does not 
differ significantly from the earlier stages. As we demonstrate 
below, for SIDM there is also a transition between two evolu- 
tionary stages, but in this case the difference is very significant. 
These two stages are best described by distinguishing between 
two limits: 

The Long Mean Free Path (Imfp) Limit. If the typical distance 
a particle travels between collisions is much longer than the 
gravitational scale height at the particle's position, it orbits the 
cluster center many times unperturbed before being scattered. 
For typical parameters, the low-density, extended halo of a ther- 
malized SIDM system will reside in the Imfp regime, while the 
core may or may not. If the core is in the Imfp regime, a parti- 
cle can escape from the core following a single encounter with 
another particle, since momentum transfer in these large-angle 
encounters is appreciable. The entire volume of the core can 
then participate in heat and mass transfer. 
The Short Mean Free Path (smfp) Limit. If the core is dense 
enough, its size may be comparable or even exceed the mean 
free path between collisions. When the mean free path is much 
smaller than the gravitational scale height, particle motion in 
the core is constrained by multiple collisions. In particular, 
heat conduction is obviously less efficient, and transfer of en- 
ergy and mass from the core to the halo is limited to a sur- 
face effect at the edge of the core. As demonstrated below, 
even if a SIDM halo is formed initially with a core in the Imfp 
regime, gravothermal evolution will cause the central density to 
increase and ultimately drive the core to a smfp state. We also 
show that in this smfp limit, mass loss from the dense core is 
severely reduced with respect to a Imfp core. 

In this paper we apply the gravothermal model to study the 
evolution of an isolated, spherical halo of SIDM. Although this 
is clearly a simphfied description omitting many of the details 
of structure formation (mergers, accretion, angular momentum 
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transfer, baryon content), our approach offers a simple tool for 
assessing the SIDM proposal, and provides some new predic- 
tions regarding the evolution and late-time gravothermal col- 
lapse of such a system. The model also provides physical in- 
sight for interpreting the results of more detailed N-body sim- 
ulations, and serves as a probe for studying details of the core 
profile not easily resolved by these simulations and the late- 
time behavior whose diagnosis is beyond their current reach. 

In § 2 we estimate relaxation timescales characterizing typ- 
ical SIDM halos and compare them with those of globular star 
clusters. In § 3 we outhne the gravothermal formahsm. The 
results of our SIDM halo evolution calculation are presented in 
Sections 4 and 5. We first show in § 4 that in the Imfp limit 
there exists a self-similar solution for the halo evolution, equiv- 
alent to the Lynden-Bell & Eggleton (1980) homology solution 
for star clusters. Then, in § 5 we perform a time-dependent 
gravothermal numerical calculation to examine the full evolu- 
tion, including the later stages where the core moves into the 
smfp regime. We briefly compare our results with previous N- 
body simulations in § 6, and present our conclusions and a dis- 
cussion in § 7. 

2. TIMESCALES 

In a globular star cluster, the thermal relaxation time due 
to multiple, small-angle, gravitational (Coulomb) encounters 
is given by (see, e.g., Lightman and Shapiro 1978 and Spitzer 

1987) 

triGC) = 



15 A(?m^nlogiO AN) 



; 5 X lO^rs 
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1/2 



-1/2 



5 pc 



3/2 



(1) 



In equation (1) A' is the number of stars, m is their mass, n is 
their number density, and v^ is the stellar velocity dispersion, 
related to the total cluster mass M and the half-mass radius, Rh 
by the virial relation v^ « \{GM/Rh). In typical virialized star 
clusters, the ratio of collision mean free path to local system 
scale height H everywhere satisfies the strong inequality 

HJgc gc ~ 261og(0.4A0 ^ ^ ' 

where td w H /\m is the local crossing or dynamical timescale. 
As a result of this inequality, these systems are only weakly col- 
lisional, satisfying the coUisionless Boltzmann (Vlasov) equa- 
tion to high approximation as they evolve quasistatically on 
relaxation timescales. For typical globular cluster parameters, 
the relaxation time is significantly shorter than the cluster age 
(« 10'*^ yrs), hence all clusters are well thermalized. Indeed the 
timescale for gravothermal core collapse to a singular state in 
such systems satisfies tcoiiiGC) ~ 330?r(GC) (Cohn 1980) and is 
comparable to typical cluster ages, so that many, initially dense 
clusters are likely to have undergone complete gravothermal 
collapse (Press and Lightman 1978; Goodman & Hut 1985; 
Spitzer 1987). 

In an SIDM halo where thermalization is due to close, large- 
angle interactions, the relaxation time is the mean time between 
single colUsions. For a cross section per unit mass a, this time 
is 

?r(SIDM) = ^ 



apva 
1.40 X 10 Vs 
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10-24gm cm-3 
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lO'^cm 



sec 
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cm^ gm 



where Pc is the central density and v^ is the central velocity dis- 
persion. The constant a is typcially of order unity; specifically, 
a = sJldj-K « 2.26 for particles interacting elastically like bil- 
lard balls (hard spheres) with a Maxwell-Boltzmann velocity 
disfi-ibution (Reif 1965, eqs. (7.10.3), (12.2.8) and (12.2.12)). 

As discussed below, for a typical SIDM halo, the ratio XjH 
is likely to be > 1 at formation throughout the halo; with evolu- 
tion, this ratio does not change much in the halo, but ultimately 
drops well below unity in the central core. At this point the 
system consists of a fluid core surrounded by a weakly colU- 
sional extended halo. By construction, the typical SIDM halo 
relaxation time is smaller than a Hubble time, enabling the sys- 
tem to thermahze. We show below that a typical halo has a 
gravofliermal collapse time tcoiiiSIDAf) w 290tr(SIDM). Ac- 
cordingly, halos forming in the present epoch are in no danger 
of collapsing, but halos formed earlier at high and moderate 
redshift may have already undergone gravothermal collapse. 

The assumption of hard sphere interactions adopted here pro- 
vides the simplest description of dark matter particle interac- 
tions and is the basis of most N-body studies to date. How- 
ever, more complex interactions, characterized by cross sec- 
tions which may be energy dependent and/or non-elastic, are 
also possible and may yield different evolutionary tracks than 
the ones we will derive. 

3. THE GRAVOTHERMAL MODEL FOR EVOLUTION OF A 
SIDM HALO 

The fundamental equations describing a spherical, viriahzed 
gravothermal fluid are mass conservation, hydrostatic equilib- 
rium, an energy flux equation and the first law of thermody- 
namics. The latter equation is time-dependent and drives the 
quasistaticevolution. Denoting M(r) as the mass enclosed by 
radius r, p{r) as the density, v(r) as the one-dimensional veloc- 
ity dispersion and L(r) as the luminosity through a sphere at r, 
these equations are (Lynden-Bell & Eggleton 1980): 
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dt 
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(5) 



(6) 



(7) 



dt 



log 



'M \P . 

We can define a pressure, p = pv^, which appears in the equa- 
tion of hydrostatic equilibrium, and a temperature ksT = mv^ 

which appears in the flux equation; here m is the particle mass. 
The last equahty in equation (7) introduces an effective entropy. 



(8) 



5 = log — 

P 



Note that the time derivatives in equation (7) are Lagrangian. 
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The detailed form of the flux equation depends on the nature 
of the heat conducton. In standard heat diffusion, where the 
mean free path between coUisions is significantly shorter than 
the system size, the flux equation is 

where is is an effective "impact parameter" of order unity and 
A is a coUisional scale for the mean free path given by 

A= — . (10) 

pa 

For a gas of hard spheres with a Maxwell-Boltzmann distri- 
bution, the mean free path is actually (. = X/Vl (Reif 1965, 
eq. (12.2.13)). The coefficient b can be calculated to good 
precision from transport theory, and has the value of b ^ 
2577/(32^6) f« 1.002 (Lifshitz & Pitaevskii 1979, chapter 1 
eq. (7.6) and problem 3). Note added in proof - a similar for- 
mualtion for the fluid limit of SIDM was presented by Gnedin 
& Ostriker (2001) who studied substructure in galaxy clusters. 

Equation (9) does not apply for a dilute gravothermal system 
where the mean free path is significantly larger than the gravi- 
tational scale height of the system H, where 



H: 



AnGp 



1/2 



(11) 



In this case, particles make several orbits between collisions, 
varying their radial position by about one scale height. If 
X^H, the mean time between coUisions, tr, is then larger than 
the dynamical time, td=H/y, and the flux equation in this Imfp 
limit can be approximated as 

L 3, H^dy^ 

(Lynden-BeU & Wood 1968; Spitzer 1987). The two forms of 
the flux equations can be combined into a single expression in 
order to treat a general gravothermal system, obtaining 

This equation correctly reduces to Equation (9) in the smfp 
limit and to Equation (12) in the Imfp limit. Substituting in 
Equation (4) and Equation (10) yields 



-^=--b V 

47rr2 2 
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dr 



(13) 



- — T = —-abva 
4nr^ 2 



AttG 



pv^ 



dr 



(14) 



This combined form of the flux equation arises from our need 
to calculate heat diffusion from the dense core to the dilute, ex- 
tended halo. While the extended halo always resides in the Imfp 
Umit, the character of the core may change. Consider the ratio 



A = (l/pa) (vV47rGp) 



-1/2 



= (4ttG) 



1/2 



Vp) 



-1/2 



(15) 



Assume that initially the core is dilute enough so that imme- 
diately following virialization it resides in the Imfp limit where 
the ratio is above unity. During gravothermal evolution the core 
parameters Vc and both increase, decreasing the ratio \/H. 
Now the core density increases much faster with time than the 
core temperature (dlog{vlit))/dlog{pcit)) « 0.1, see below). 
Hence, a nonrelativistic core will enter the smfp regime well 
before gravothermal contraction drives it to a relativistic state 
and dynamical collapse. Therefore, late in the evolution of a 
SIDM halo, its core will be in the smfp regime, while the ex- 
tended halo will remain in the knfp regime and some transition 
region must exist between the two. 



4. SELF-SIMILAR SOLUTION FOR EVOLUTION IN THE LMFP 

LIMIT 

The equations (4, 5, 7 and 14) describe the evolution of 
a spherical, isolated halo of SIDM. These equations define a 
time-dependent diffusion problem, which we solve numerically 
in the next section. However, it is useful and instructive to ex- 
amine first a self-similar solution which can be found for these 
equations in the Imfp limit, following the original derivation of 
Lynden-Bell & Eggleton (1980) for globular clusters. 

If the entire halo is sufficiently dilute so that H every- 
where, then flCT^ can be neglected with respect to 47rG(/3V^)"' 
in the square brackets in equation (14). This Imfp regime ad- 
mits a self-similar solution which can be obtained by separating 
variables into time and space components as follows: 

r = rc{t)r^ ,M = Mc(t)M^ , v = Vc(Ov* , p = pc(t)p* ,L=Lc(t)L^ . 

(16) 

The time-dependent functions give the values of the core pa- 
rameters at time t. The dimensionless spatial profiles denoted 
by the asterik (*) are related by a set of ordinary differential 
equations, found by substituting Equation (16) into the four 
equations (4-7): 



(17) 
(18) 
(19) 

C1-C2 7; • • (20) 



dr^f 


= rip* 


dip*yl) 


M^pt, 


dr^ 


r^ 


d{vl) 




dr^ 


dp*^l 


dL^, 
dr^ 


2 2 



dlogjvl/p^,) 
JlogM* 



The constants ci and ci are two eigenvalues of the equations. 
The system of four equations plus two eigenvalues is uniquely 
determined by specifying six boundary conditions. Four of 
these conditions apply to the origin, where 



/9*(0) = 1, v*(0) = 1, M*(0) = aiidL*(0) = , 



(21) 



and the other two conditions arise from the requirement that 
the extended halo remains practically static over the evolution 
of the core, so that it may be treated as infinite. The mass 
and energy of the extended halo are then also infinite, elimi- 
nating any scale to the problem and ensuring the existence of 
a self-similar solution. The consequence of this condition on 
the extended halo is that its density profile must tend asymp- 
totically toward a power law of the form /9*(r* 00) oc r~" 
(Spitzer 1987). Indeed the condition of self- similarity automat- 
ically yields the relation a = 6(c\ -C2)/(2c\-C2) (Lynden-Bell 
& Eggleton 1980). The asymptotic outer boundary conditions 
then arise from equations (17, 19) and are 



Z.*(r* ^ 00) = (a-2)/9*v^r* 



v*(?-* 



00) : 



r*(2a-2) 

(22) 

The above derivation is almost identical to the derivation 
of the self-similar solution for star clusters by Lynden-Bell & 
Eggleton (1980). The sole difference is the form of the flux 
equation (19), since for SIDM tr oc (pv)"\ whereas for globular 
clusters tr ocv^/p. 
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not taken into account). This relation is 



Fig. 1 . — The self-.similar solution for SIDM in the Imfp limit. Mass(M»), 
density (p* ), ID velocity dispersion (v* ), and luminosity (L* ), are plotted as 
functions of radius (r«), all in nondimensional units. 

We solved the coupled set of ODEs for the nondimensional 
profiles by relaxation on a finite mesh (Eggleton 1971), and the 
solution is shown in Figure 1 . A key feature of the solution is 
the existence of a core-halo structure: the finite core is nearly 
isothermal and homogeneous, while the infinite halo is charac- 
terized by declining power-law profiles for the density, velocity 
and luminosity. Note that the luminosity peaks at the outer edge 
of the core (r* — 1), reflecting the shorter evolutionary (cooling) 
timescale (~ M^/ (Lr)) in this region. 

A fundamental feature of the self-similar solution is that the 
eigenvalues ci and C2 completely define the time-dependent 
evolution in terms of pc(0),Vc(0). We can define an initial core 
radius rdt = 0) and mass, and Mdt = 0) and follow their time 
dependence as well. Delineating the precise core boundary for 
the self-similar density profile plotted in figure 1 is somewhat 
arbitrary. We will employ the King radius (King 1966), where 

= 3H{pc,Vc) and then define Mc = M{rc). Typically, the den- 
sity at the King radius is about 1/3 of the central density. 

The time-dependent evolution can be parameterized by a sin- 
gle variable (, which is defined through 



dE, _ E, dM, 
dt Mr dt 



(23) 



(Spitzer 1987) where Ec is the total energy of the core. Simple 
dimensional arguments lead to the following relations between 
core mass and other core quantities: 

^logr.(0_2-C, ^^^^=-(1-0, ^^"^^ =-(5-30: 



dlogMcit) 



dlogMcit) 



dlogMcit) 



and C is directly related to the exponent a according to 

5-2a 



(24) 



(25) 



The self-similar solution yields the time dependence of the core 
quantities in terms of the "coUapse time", tcoii(0), the time (from 
t = 0) when the core mass and radius vanish while the density 
and velocity dispersion go to infinity (relativistic instability is 



Mc(t) 
MM 



1-- 



tcolliO), 



where 



and 



tcoim=e-trio), 



11-7C 



(26) 
(27) 
(28) 



The mass evaporation rate parameter, ^g, is defined through the 
relation 

dMrit)_ MS) 



dt 



tr(t) ' 



and is given by 



6 = 2'^2b , 



(30) 



for the Imfp flux equation (12). Similar time-dependent expres- 
sions for the evolution of r^, Vc, and pc follow from equations 
(24) and (26). The sole difference between this solution and 
the star cluster case arises from the different functional form of 
the relaxation time: for Coulomb scattering 9 = 2/(7-30- In 
Table 1 we list aU the relevant parameters for the self-similar 
SIDM solution and compare them with the star cluster solution 
(Lynden-BeU & Eggleton 1980; Spitzer 1987). 

Most of the parameters for the SIDM evolution problem have 
values which are very similar to the star cluster case. This is 
mostly due to the fact that the allowed range of the parameters 
is quite small. SpecificaUy, the exponent a must reside between 
a = 2.0 (^ ( = 1), which corresponds to an isothermal profile, 
and a = 2.5 (^ C = 0) which corresponds to a constant core 
energy (see Lynden-Bell & Eggleton 1980 for a discussion of 
this allowed range). In star clusters, transfer of particles from 
the core to the extended halo is a diffusive process, so the av- 
erage change in specific energy should indeed be smaller than 
in SIDM systems, where a particle can be ejected from the core 
by a single colhsion. Hence, the star cluster values set an up- 
per (lower) limit on the value of a (Q for the SIDM case, as 
is borne out by the actual solution. An immediate consequence 
is that the density of dark matter in an isolated, extended halo 
of an SIDM system should maintain the profile p(r) oc r"^-^^. 
Such a power law is not inconsistent with halo profiles previ- 
ously found with N-body simulations for SIDM. In the same 
way, the predictions of the self-similar gravothermal model for 
star clusters have been confirmed by detailed integrations of the 
Fokker- Planck equation (Henon 1961; Spitzer and Hart 1971; 
Cohn 1979; Marchant & Shapiro 1980). 

The result for ( corresponds to the core mass evolving with 
velocity according to the relation dlogMc/dlogv^, » -4.27. 
The core of a Imfp SIDM halo must lose more than 99.99% of 
its mass in order to increase its central temperature by one order 
of magnitude. Assuming that the core is initially nonrelativis- 
tic, it must increase its central temperature by a factor of, say, 
10"*- lO*" to reach the relativistic instability, in which case only 
a tiny fraction of the initial core mass will be available for the 
black hole that arises from the final catastrophic collapse. We 
show in § 5 below how this conclusion is drastically modified 
once the core enters the smfp regime. 

One key issue concerning galaxy formation is how the col- 
lapse time compares with typical ages of galaxies. According 
to Table 1, a Imfp SIDM system has a lifetime of about 290?r 
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Table 1 

Parameters Dehning Self-similar Gravothermal Evolution for SIDM and Globular Clusters (GC). 





a 


c 


C\ 




Q 


6 


(/loi5(M,.) 


dlogCft.) 
rflog(M,) 


tcoU 
tr 


SIDM 


2.190 


0.7655 


1.903x10"^' 


8.092x10"^ 


0.3545 


1.21x10-3 


-0.2345 


-2.704 


-291 


GC 


2.208 


0.7382 


2.322x10-3 


9.704x10-4 


0.4179 


1.31x10-3 


-0.2618 


-2.785 


~319 



before undergoing secular core coUapse. This ratio tcoii/tr is 
somewhat smaller than in the star cluster case, because of the 
greater efficiency of energy transfer in clusters as reflected in 
the value of b. The cumulative effect of multiple gravitational 
scatterings in the star cluster case gives b « 0.675 (Spitzer 
1987), whereas for hard sphere colUsions in SIDM systems, we 
have b=l .002. The typical density of dark matter inferred for 
flat-density cores is — 0.02 Mq pc-^ (— 1 .4 x 10"^"* gm cm"^ ) 
(Firmani et al. 2001), so for v^ < 10^ cm s^' the collapse time 
exceeds ten Hubble times (for cr ~ 1 cm^ gm-^ ). Consequently, 
SIDM is consistent with the existence of relaxed, flat-core den- 
sity profiles in present day halos. Furthermore, the functional 
forms of Equations (24) and (26) show that throughout most 
of the evolution, the core hardly changes: at / = 0.5feo;/(0) the 
central density wiU have increased only by a factor of — 2. Cur- 
rently existing flat-density cores should exhibit central densities 
similar to those they acquired at virialization. 

5. TIME-DEPENDENT EVOLUTION OF A SIDM HALO 

We now return to a full time-dependent, numerical solu- 
tion of the evolution equations of § 3, without assuming self- 
similarity or restricting ourselves to the Imfp regime. The nu- 
merical approach is not only useful for confirming the self- 
similar solution, but it is essential for studying the evolution of 
the core in the smfp limit, which it must eventually reach. Late 
in its life (shortly before tcn), the core becomes sufficiently 
dense and hot that A/7/ falls below unity and the evolution devi- 
ates from the self-similar solution. This final stage will be dom- 
inated by the thermal evolution of a genuine fluid core, which 
cools, contracts and ultimately becomes dynamically unstable 
to relativistic collapse to a black hole. 

It is useful to define a new set of dimensionless variables us- 
ing fiducial mass and length scales, Mq and Rq. The interaction 
cross section can then be written as 



47r/?g 

(T = (TCTo, With (To = 

Mo 

The natural scales for the other variables are 



Mo 



, vo = 



1/2 



GMl 
Roto 



We define the characteristic time scale, to as 

Ro _ tr.O 
ba ' 



to- 



ol) 



(32) 



(33) 



abvQcr 

where f^.o is the relaxation time for po, vo and ao. This choice 
allows to us to recast the evolution equation in a dimensionless 
form: 

dM 



df 
df 



Mp 



(34) 



(35) 



9(vVp) 3 1 d 
dt ~ Ipdf 



^2 1 
pv^ 



df ■ 



(36) 



In deriving equation (36) we combined the flux equation (14) 
and the first law of thermodynamics (7) to eliminate the lu- 
minosity and arrive at a diffusion equation. The tilde quanti- 
ties represent the dimensionless variables, i.e., x = x/xq. Note 
that unlike the self-similar solution, these generalized equations 
are not scale-free due to the presence of the aa^ term in equa- 
tion (36). 

By hypothesis, we assume that a SIDM halo is initially dilute 
enough to be in the Imfp hmit, so we can use the the self-similar 
profile as the initial condition. The natural choices for scaling 
are then vo = Ve(f = 0) and po = pc(t = 0), so that Vc(f = 0) = 1 ,and 
Pcit = 0)= 1. Hereafter we drop all tildes but report results in 
terms of these dimensionless variables. 

We solved this time-dependent diffusion problem nu- 
merically by taking alternative diffusion and hydrostatic- 
equilibrium steps. Time dependence is followed through equa- 
tion (36), and after each successful time step the profile is ad- 
justed to maintain hydrostatic equiUbrium. Since the dynami- 
cal time scale of the system is much shorter than the thermal 
time, hydrostatic equilibrium is solved by keeping the entropy 
of each Lagrangian zone fixed, while its position and other ther- 
modynannic properties are modified. The effect of the SIDM 
collisions is examined by varying the value of the combina- 
tion fl(T^. In principle such a (Newtonian) calculation can be 
carried out all the way to collapse, but in practice maintain- 
ing numerical accuracy requires constantly refining spatial and 
temporal resolution, which ultimately hmits the dynamic range 
of the calculation. Reasonable accuracy can be maintained un- 
til to v^/v^(f = 0) » a few tens, at which point theb integration 
time step become excessively small. We are able to asses the 
approximate state of the core at the onset of the relativistic in- 
stability by extrapolating the results of our late integrations to 
the relativistic epoch. 



5.1. Confirmation of the Self-Similar Solution 

Our hydrodynamic integrations confirm the self-similar so- 
lution for the Imfp limit found in § 4. By setting aa^ = in 
equation (36), the entire halo is forced to reside in this limit. 
We show in Figure 2 snapshots of the density profile p{r) for 
aa^ = 0. The initial profile is marked by the thicker line, and 
in the other profiles larger central densities correspond to later 
times. The self-similar structure is clearly maintained through- 
out the evolution, as the core contracts in size and mass, while 
the extended halo remains largely unchanged. Numerically we 
find that the density in the extended halo falls as p cx ^--^ 
corresponding to a value of C — 0.766. The quantitative fit to 
the self-similar solution is very good, as is the fit for the func- 
tion pc{t) and other quantities (see below). 



SIDM AND GRAVOTHERMAL CATASTROPHE 



7 




r/r,(0) 

Fig. 2. — Snapshots of the density profile of a Imfp halo of SIDM at 
selected evolutionary times. The thick line denotes the profile at ? = 0; 
larger central densities correspond to later times. Profiles are drawn at f/to = 
0.0, 264.0, 286.5, 288.52, 288.70, 288.7178 and 288.7180. 

5.2. Evolution to a Short Mean Free Path (smfp) Core 

We now consider our results for aSP' > 0, whereby the core 
eventually evolves into the smfp regime. Figure 3 shows snap- 
shots of the density profile for a calculation with acP- = 0.01; 
the profiles are labeled according to the value of X/H (in units 
of a^^^^) at the center. At early times, when X/H > 1 every- 
where, the profiles are identical to those in figure 2, but as the 
core density and temperature increase, the density profile grad- 
ually deviates from the self similar solution. Some difference 
is already observable when X/H < 0.1, and when the core has 
A/i/ < 10"^ the density profile has obviously taken a new form. 
In fact, what was originally the core during the Imfp stage of 
the evolution has fragmented into two components: a dense in- 
ner core which continues to evolve, and an outer region with a 
lower, roughly uniform density, which connects to the extended 
halo. This outer core is also almost static while the inner core 
continues to contract with time. 

The existence of the double core structure is even more 
pronounced when examining the evolution as a function of 
the mass coordinate. In Figure 4 we show snapshots of the 
density and temperature profiles for aa^ = 0.01 (the times of 
the snapshots are not necessarily identical between the two 
frames). The inner core develops a fairly sharp "edge", which 
becomes more distinct as the evolution progresses. The na- 
ture of this inner core can be understood from the (X/H)(m) 
and L(m) profiles at the late stages of the evolution. These are 
shown in Figure 5 (along with v^(m)) for the last profile of 
Figure 4. The inner core forms in the region where H ^ X, 
and heat conduction is dominated by classical diffusion with 
a short mean free path. Its edge corresponds to a transition 
region where // ~ A and where the luminosity peaks. The in- 
ner core resembles a fluid star composed of a nonrelativistic, 
nearly isothermal and homogeneous gas undergoing thermal 
cooling, surrounded by a thin, exponential envelope. The outer 
part of the core is the remnant of the original core in the Imfp 
limit. This outer core then coimects smoothly to the extended 
halo, and barely evolves with time as the inner core contracts. 




r/r,(0) 

Fig. 3. — Snapshots of the density profile of a SIDM halo with 
aa^ = 0.01 at selected evolutionary times. The thick line denotes the pro- 
file at t = 0; larger central densities correspond to later times. Profiles 
are drawn at (l/lo,X/Ha^/^) = (0.0,10.0), (264.0,2.8), (287.0,8.2 x lO'^), 
(289.5,2.1 X 10-1), (289.90,6.0 x lO'^), (290.12,1.3 x lO'^), (290.30,7.0 x 
10-3), (290.49,3.0 X 10-^), (290.63,2.1 x 10-^). 

An important aspect of the composite core is that the outer 
region holds most of the mass, and therefore acts as a heat sink 
to the inner region, which continues to evolve as its outer lay- 
ers cool. The iimer core loses mass as it continues to contract. 
Gravothermal collapse is slowed down with respect to a hnfp 
regime, where particles orbit outside the core for many scale 
heights; in a smfp core, mass and heat loss are now confined 
to the surface where A « //. A slight temperature inversion 
is created in the outer core, which corresponds to a region of 
negative (although very small in magnitude) luminosity. The 
timescale for heat conduction in the outer core is much longer, 
however, than the thermal time of the inner core, so the details 
of the outer core temperature profile are unimportant for the 
gravothermal evolution. 

The value of the cross section, a, determines the evolution- 
ary track of the SIDM halo as its core contracts. The deviation 
from the self- similar solution wiU begin earlier (at lower cen- 
tral densities and temperatures) for larger a; see equation (14). 
A larger cross section also reduces the efficiency of heat con- 
duction, since the conductivity scales as . These effects are 
demonstrated in Figure 6, which presents the time dependence 
of the central density of a SIDM halo for different values of 
a&^. In all cases the initial conditions are the self-similar so- 
lution. In the extreme Imfp limit (aa^ = 0), the collapse time 
appears to be ~ 289fr(0), in good agreement with the predic- 
tion of the self-similar solution (assuming b ~ 1). With a finite 
aa^ the late-time evolution is slower (requires more relaxation 
times), as expected, and the delay of the late-time evolution is 
more pronounced for larger values of the cross section. 

A finite value of the cross section introduces a scale into the 
problem which precludes a simple, self- similar solution for 
the evolution of the double core and extended halo structure. 
For example, we show in Figure 7 the relation between central 
density and central temperatures for calculations with different 
aa^. The details of the evolution depend on the value of a: 
there is no universal pdy^) relation. We note in passing that 
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even if the entire system (core and extended halo) were in the 
smfp limit, a self similar solution with a power law density pro- 
file cannot exist: the thermal time scale in the core would be 
longer than that in the halo, causing the system to be unstable 
(Hachisu et al. 1978). 
10^ 





10 

m/M,(0) 



270 275 280 285 290 295 300 305 
t/t„ 

Fig. 6. — Central density as a function of the time for a Imfp SIDM halo 
(aa^ = 0) and for SIDM halos with a positive value of aa^ . The main figure 
shows the late time evolution, when central density increases rapidly; the inset 
shows the entire evolution from t = 0. 



Fig. 4. — Snapshots of the (a) density and (b) temperature profiles of a SIDM 
halo with aa^ = 0.01 as functions of the mass at selected evolutionary times. 
The thick fine corresponds to the initial conditions. Note the formation of an 
inner core with a relatively sharp surface. 



10 



(a) 




10' 



Fig . 5 . — The (a) temperature, (b) luminosity and (c) relative mean free path 
\/H profiles of a SIDM halo with aa^ = 0.01 at f = (thick line) and at the 
termination of the calculation, t/tr(0) ~ 290.713 (thin line). 




vVv'e(O) 

Fig. 7. — Central density vs. central temperature for a SIDM halo. Curves 
are labled with the value of air^ used in the simulation The Imfp self-similar 
solution is plotted for aa^ = 0. 

The most significant effect caused by the the transition of the 
irmer core into a strongly coUisional fluid is a critical reduction 
of the rate of mass loss from the core. We demonstrate this in 
Figure 8, which shows the core mass as a function of the cen- 
tral temperature for the different cases. In the Imfp phase, the 
core mass follows a power law with <ilogMc/<ilog(Vp) ~ -4.27, 
in agreement with a value of C = 0.766 (see eq. (24)). When 
0(7^ > 0, the mass loss rate is appreciably reduced for the smfp 
inner core: instead of a steep power law, the core mass nearly 
saturates while the central temperature increases. As a con- 
sequence, the core mass at the onset of relativistic instability 
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should be much larger than would be expected from the self- 
similar solution in the Imfp limit. While strict self- similarity is 
broken once the inner core enters the smfp regime, there does 
seem to exist some power-law correspondance for the evolu- 
tionary tracks of the central values of some of the physical pa- 
rameters. Empirically, we find that the function M^ja^'^ (in 
units of a""^) tends to converge to a single track as a function 
of central temperature; this convergence is also shown in fig- 
ure 8. 



;>\ \ 

- 10"' 










\ '"-SV—i^^o 
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10^ 




~~~~ 1/30 




ao'=10^ 











Ve'/V\(0) 

Fig. 8. — The core mass vs. the central temperature for a SIDM halo (solid 
lines). Curves are labled with the value of au^ used in the simulation. The 
Imfp self-similar solution is plotted for acf^ = 0. Also shown are curves nor- 
malized as Mcliac^'f -^ for the cases with acP' > 0. (dashed lines). 

5.3. Late-Time Evolution of a Smfp Core 

The consequences of relativistic instabihty in the inner core 
depend on the colUsional nature of the gas at the onset of in- 
stability. For example, when a fluid star is unstable to radial 
collapse on a dynamical timescale, it undergoe catastrophic col- 
lapse in which the entire star (core + halo) forms a black hole 
(see, eg. Shapiro & Teukolsky 1980). By contrast, when a col- 
lisionless cluster of particles with a core-halo structure is un- 
stable to radial collapse, it is essentially the core alone and its 
immediate surroundings that undergoes catastrophic collapse 
(Shapiro & Teukolsky 1986,1992). Particles in the halo out- 
side the core constitute the bulk of the matter and continue to 
orbit the central black hole, forming a new, dynamically sta- 
ble, equilibrium system consisting of a central black hole and 
a massive, extended halo of orbiting particles. When a SIDM 
cluster becomes sufficiently relativistic in the inner core it too 
must be subject to relativistic instabihty to collapse, as in the 
two extreme opposite cases described above. At the onset of 
collapse, the inner core acts as a fluid system (A/// <C 1) while 
the exterior region behaves as a coUisionless gas (X/H ^ 1), 
at least on dynamical timescales. This bifurcation is evident in 
Figure 5. The key consequence is that we expect the entire in- 
ner core of a SIDM cluster to undergo collapse to a black hole, 
while the region outside the inner core relaxes to a dynamically 
stable equihbrium system of particles that continue to orbit the 
central hole. 



Thus the relevant mass to consider for dynamical collapse 
and black hole formation is the inner core, or, more specifically 
the region in which the matter satisfies X/H <l and behaves 
as a fluid. The matter inside will all collapse and the matter 
outside this region will reach a stable equilibrium slate outside 
the central black hole, at least on dynamical timescales. Subse- 
quent interactions between particles in this ambient region wiU 
ultimately fuel the central hole, causing it to grow on relaxation 
timescales (Ostriker 2000). The same effect drives the secu- 
lar growth of black holes immersed in ambient star clusters, a 
problem that has been well studied by detailed Fokker-Planck 
calculations (Bahcall & Wolf 1976, Frank & Rees 1976, Light- 
man and Shapiro 1977, Cohn & Kulsrud 1978, Marchant & 
Shapiro 1980; see Shapiro 1985 for a review and references). It 
is the catastrophic collapse of the inner (fluid) core of a SIDM 
cluster which naturally provides the initial seed black hole. 

The Mc vs. v^ relation for the late-time inner core evolution 
determines the mass of the core at the time of relativistic insta- 
bility. Here we infer that relation by a combination of physical 
reasoning and extraplolaton from our numerical integrations for 
earlier epochs into the relativistic regime. Because energy and 
mass transfer from the smfp inner core are Umited to surface 
phenomenon, the relative rates of these processes must remain 
significantly reduced with respect to the Imfp core. In equa- 
tion (23), the typical absolute value of C, during the late time 
evolution of the inner core must be close to zero, since only a 
small fraction of the inner core energy is transferred to the outer 
core. If 1^1 ~ then we have (ilogv|:/(:/logM(. « -1: hence, at 
late times the inner core will lose about one order of magnitude 
in mass for every order of magnitude it gains in central central 
temperature (note that C, cannot be identically zero, as this value 
corresponds to no evolution whatsover). 

We can calibrate this prediction for the late-time Mc vs. v^ 
relation by calculating the value of Q numerically using Equa- 
tion (24). In Figure 9 we show the corresponding estimates of 
C = d\og\l/d\ogMc+ 1, using the results shown in Figure 8. 
The deviation of the time-dependent value of C, from the con- 
stant value of 0.766 found for the Imfp limit is clear and reflects 
the formation of the smfp inner core. During its formation, the 
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Fig. 9. — The exponent C derived numerically from the relation between 
core mass and the central temperature for a SIDM halo, following Figure 8. 
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mass of inner core is almost constant while its density and tem- 
perature profiles are adjusted, and so (, decreases to a minimum 
(negative) value. However, once the inner core is sufficiently 
dense, mass is continuously lost from its surface as outer layers 
cool and expand to join the outer core. The resulting value of ( 
starts to increase again, and is still rising when the calculations 
are terminated. 

We find that due to the sensitivity of the core mass to the 
details of the density profile at the edge of the inner core, it 
is preferable to estimate the asymptotic trend of C by a more 
robust definition of the core mass. In the context of the dou- 
ble core structure where the inner region is coUisional (fluid), 
the transition point where X/H = 1 provides a more appropriate 
definition for the "edge". In Figure 10 we again show the nu- 
merical estimate of (, this time when the core mass corresponds 
to the radial position where X/H = 1 (which only exisits after 
a smfp core is formed). Again we see that for all cases C, is in- 
creasing from a transient stage of C ~ -1, and seems to settle 
asymptoticaUy in the range -0.2 < ^ < -0.1. We expect that 
this trend is indeed representative of the typical value of ( in 
the case of a smfp core at high central temperatures. Our pre- 
diction is, therefore, that once a smfp inner core is completely 
formed, the final stages of its evolution will be characterized by 
Jlog(v^)/Jlog(Mc) « -1.1 ^ -1.2. The iimer core mass wiU 
therefore decrease by only about one order of magnitude for 
every order of magnitude that its central temperature increases. 

Finally, our inability to integrate the full system of equations 
"forever" does not prevent us from estimating the the collapse 
time, tcoii- In the Imfp regime, the collapse time is always 
~ 290 times the instantaneous relaxation time. At the time at 
which the numerical simulations are halted, the smfp inner core 
relaxation time is typically 10^*^ times shorter than the initial 
core relaxation time. Consequently, the total time to collapse 
from f = is dominated by the early evolution of the Imfp core 
and extended halo system. The additional time to collapse once 
the inner core enters the smfp regime is considerably shorter. 



0.0 




v>A0) 



Fig. 10. — The exponent ( derived numerically from the relation between 
core mass and the central temperature for a SIDM halo, when the core mass is 
defined by the position where \/H = 1 . The asymptotic value of = —0. 175 at 
late times is indicated approximately by a horizontal dashed line. 



6. COMPARISON WITH N-BODY SIMULATIONS 

Our analysis of SIDM evolution via a gravothermal fluid ap- 
proach is an alternative to N-body simulations. In general, 
mean-field N-body codes have the advantage that they avoid 
having to adopt an approximate, fluid-like set of equations 
to describe SIDM. However, the spatial resolution of N-body 
codes is severly hmited by the total number of particles they 
are able to follow. Large numbers of particles, required for high 
resolution, can be tracked only for restricted integration times, 
and this hmitation prohibits N-body simulations from following 
the evolution of SIDM halos for many relaxation timescales. 

Previous N-body calculations focussing on the same problem 
tackled here (Burkert 2000; Kochanek & White 2000; Yoshida 
et al. 2000; Dave et al. 2001) also found that introducing col- 
Usions in the dark matter interactions leads to the formation of 
a flatter core, as opposed to the CDM prediction of a cuspy 
core (with p(r) cx r"", and a = -1 in Navarro et al. 1997 or 
a ~ -1.5 in Moore et al. 1999). Previous simulations of time- 
dependent evolution also confirm the general trend of core con- 
traction while the extended halo remains roughly static. They 
also show that during the evolution the density and temperature 
of the core both increase (see, e.g. in Fig. 1 of Burkert 2000, 
Figs. 1 and 2 in Kochanek & White 2000 and Fig. 3 in Yoshida 
etal. 2000). 

Our study is compared most naturally to the simulations of 
Burkert (2000) and Kochanek & White (2000), who also inves- 
tigated an isolated halo. Kochanek & White (2000) surveyed 
different values for the SIDM cross section, and discovered that 
a larger cross section induces a greater effect on the evolution 
- see their Fig. 2b, which resembles our Figure 6. However, 
both works assumed as initial conditions a virialized Hernquist 
(1990) profile with a cuspy core, so that initially, part of the 
halo resided in the smfp regime: the initial conditions of Burk- 
ert (2000) correspond to a ratio X/H « 0.1 at the outer edge 
of the core and those of Kochanek & White (2000) correspond 
to the range of 0.1-3.0. Such a system has a lifetime of only a 
few relaxation times (Quinlan 1996), which is the trend found 
by Kochanek & White (2000), who claimed that such a short 
lifetime is an argument against the validity of SIDM. However, 
such an initial configuration appears to be inconsistent with the 
SIDM hypothesis, at least for large collision cross sections. In 
this case, a cluster will not virialize with a cuspy halo as in the 
case of CDM, since coUisions would modify halo profiles. The 
formation of smooth cores following virialization is certainly 
indicated in the results of Yoshida et al. (2000) and Dave et al. 
(2001), whose calculations include halo formation by violent 
relaxation. Our results may offer a better description of the 
evolution of a SIDM halo: the lifetime of a present day halo is 
significantly longer than a Hubble time. 

Limited spatial resolution did not allow Burkert (2000) and 
Kochanek & White (2000) to follow the evolution of the core 
to large densities. In particular, the simulations they reported 
were terminated when the central density increased by less than 
a factor of 10 with respect to its value when the flat core formed. 
Thus, they could not observe the late-time evolutionary effects 
which appear at much larger densities, including the formation 
of the smfp inner core or its reduced rate of mass loss as the 
core continues to contract. We hope that N-body simulations 
or other Boltzmann-solvers with improved spatial resolution in 
the core and extended integration timescales might be able to 
confirm our findings on these issues in the future. 
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We cannot easily compare out results with those of Yoshida 
et al. (2000) and Dave et al. (2001), who examined the signif- 
icance of a SIDM model on structure formation by including 
a Monte-Carlo module in a cosmological tree code. Their cal- 
culations are three dimensional and include the effects of ac- 
cretion and mergers, which are not accounted for in our spher- 
ical study. Dave et al. (2001) claim that no halo ever devel- 
ops an isothermal core since hot material is continuously ac- 
creted onto the halo, therefore preventing efficient heat trans- 
fer from the core to the halo. Although their results cannot 
truly gauge the properties of the core due to insufficient reso- 
lution, it is clear that significant ongoing accretion will modify 
gravothermal evolution which we found above. However, the 
time-dependent results of Yoshida et al. (2000) seem to be con- 
sistent with gravothermal evolution, with the exception of the 
destabihzing effect of major mergers. We plan to extend our 
method to incorporate an effective accretion and merger algo- 
rithm in the future. 



7. CONCLUSIONS AND DISCUSSION 

The purpose of this work was to construct a simple the 
gravothermal fluid model to study the evolution of an isolated, 
spherical, halo of self-interacting dark matter (SIDM). For 
SIDM every close collision can produce a large-angle scatter- 
ing, a process which can drive thermal relaxation and gravother- 
mal core collapse. By contrast, CDM halos are collisionless and 
frozen in time following their initial virialization. For a typical 
range of SIDM cross sections and central densities of dark mat- 
ter cores at the present epoch, the relaxation time is of order 
lO' years. Accordingly, a typical halo has had sufficient time 
to thermalize and acquire a gravothermal profile consisting of a 
flat core surrounded by an extended halo. In this sense SIDM 
is similar to globular star clusters, for which the gravothermal 
model was originally developed. However, in globular clusters, 
relaxation is driven by the cumulative effect of many distant, 
small-angle Coulomb encounters. 

In the gravothermal evolution of a SIDM halo, we can distin- 
guish between a Imfp limit, where the gravitational scale height 
is much smaller than the collision mean free path, and the smfp 
limit, where the situation is reversed. Heat conduction and mass 
transfer between the core and the extended halo differ in form 
and magnitude in the two hmits. 

A semianalytic self- similar solution exists in the Imfp regime 
and can be determined following the approach of Lynden-Bell 
& Eggleton (1980) for globular star clusters. The general 
case requires numerical integration of the full set of quasistatic 
gravothermal equations. We showed that the core of a Imfp halo 
undergoes gravothermal collapse in a time tcoii = 290?^ , where 
tr is the instantaneous relaxation time. In Newtonian theory the 
core collapses to infinite density and zero mass, but in reality it 
must collapse slightly earlier to a black hole with a finite mass 
due to the dynamical instability which sets in when particle ve- 
locities become relativistic. In any case, the collapse time of 
present day halos is significantly longer than the Hubble time. 
Accordingly, the existence of halos with flat cores at the cur- 
rent epoch is entirely consistent with the SIDM hypothesis and 
gravothermal evolution. 

In the standard scenario of the gravitational collapse of an 
initially dilute, overdense fluctuation of dark matter, a nascent 
SIDM halo will virialize by violent relaxation on a dynanoical 
(collapse) time scale. The system will subsequently thermalize 



by collisions on a relaxation time scale (Eq. (4)). For typical 
parameters, both the core and extended halo will be in the Imfp 
regime following the initial relaxation. The model then predicts 
that the extended halo, if isolated, will acquire a density profile 
with p(r) oc r"^ '^; this power-law should remain constant while 
the core evolves. While a power-law of this magnitude is not 
easily distinguishable from the range considered plausible for 
CDM, p oc r-2-r-3, (Navarro et al. 1997; Moore et al. 1998; 
Kravtsov et al 1998) it is consistent with the results of N-body 
simulations of SIDM. 

While the extended halo remains almost static the core con- 
tracts and eventually enters the smfp regime where it behaves 
as a fluid. Subsequent evolution drives the core into two com- 
ponents, a dense, smfp inner core and a more dilute Imfp outer 
core. Sharp temperature and density gradients develop at the in- 
terface between the two components where A « //. Gravother- 
mal evolution of the inner core is slowed down (i.e., requires 
more relaxation times) with respect to the evolution of the orig- 
inal Imfp core. 

In the Imfp limit we find that c/logMe/c/logv^ = -4.27: the 
core mass must decrease by more than four orders of magnitude 
for every order of magnitude the central temperature increases 
(very much like the case for globular clusters). If at formation 
the core is initially nonrelativistic so that (v^/ c)^ ~ 10"* - lO"'*, 
the mass of the core at the onset of relativistic instability would 
clearly be negUgible, were the evolution to persist in the Imfp 
regime. However, the transition into the smfp state drasti- 
cally changes the evolution, so that c/logMc/iilogv^ w -0.85. 
Hence, A SIDM halo retains an appreciable inner core mass as 
it evolves towards the relativistic instability. 

Previous N-body simulations of a SIDM halo show that the 
system develops a flat core plus extended halo structure. Stud- 
ies of an isolated halo (Burkert 2000; Kochanek & White 2000) 
also find that the core evolves by contracting in size while in- 
creasing its density and temperature, as we expect based on our 
gravothermal model. The gravothermal model allows us to fol- 
low the evolution of the core to much greater central densities 
than the N-body codes, and we have able to identify the appear- 
ance of a double core and reduced mass loss from the contract- 
ing inner core. We hope that future N-body simulations with 
improved resolution in the core region might be able to confirm 
this interesting behavior. 

The most significant cosmological implication of our results 
is that an isolated SIDM halo at the present epoch with a central 
density and velocity dispersion similar to those inferred from 
observations (^ 0.02 Mq pc"^ and < 10^ cm s"' respectively) 
should have a relaxed gravothermal profile, including a flat core 
and an extended halo at the current epoch. Its relaxation time is 
less than one tenth of Hubble time, so that it had sufficient time 
to thermalize, while the lifetime of its core should exceed ten 
Hubble times, so the core should be far from collapse. Hence, 
a gravothermal SIDM halo is consistent with the observational 
inference of a flat core in dark matter halos. While our ap- 
proach does not include accretion, mergers or angular momen- 
tum, these processes can be incorporated in more complicated 
gravothermal models. 

One natural question arises: did SIDM halos formed at high 
or moderate redshift earUer in the universe have core lifetimes 
that were shorter than the Hubble time? At high red shift the 
mean density of the universe was higher, and so would be the 
typical densities of halos which form, implying shorter relax- 
ation and gravothermal collapse times. A scenario whereby the 
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gravothermal catastrophe in a SIDM halo leads to the forma- 
tion of a black hole at high red shift is intriguing as a general 
mechanism for producing supermassive black holes in galaxies 
and quasars. We address this possibility in a separate work. 
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